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1.  INTRODUCTION 

Ford  Aerospace  and  Communications  Corporation  (FACC),  Western 
Development  Laboratories  Division  (WDL)  is  conducting  a  comprehensive  program  of 
research,  sponsored  in  part  by  the  subject  contract,  on  mathematical  modeling.  The 
scientific  objectives  of  this  research  program  are  to  develop  improved  probabilistic 
and  statistical  methods  and  to  demonstrate  their  power  and  relevance  through 
application  to  specific  problems  in  stochastic  networks  and  system  reliability. 


2.  RESEARCH  RESULTS  ON  PROBABILISTIC 
METHODS 

Research  efforts  to  date  of  the  probabilistic  methods  program  at  Ford  Aerospace 
are  described  in  a  series  of  ten  technical  papers.  The  fifth  through  tenth  papers 
listed  below  report  on  research  sponsored  in  part  by  the  subject  AFSC 
Contract. 

^ .  "On  Failure  Modeling,"  A.  Lemoine  and  M,  Wenocur. 

2.  "A  Stochastic  Network  Formulation  for  Complex  Sequential 
Processes,"  A.  Lemoine. 

3.  "On  Shot  Noise  and  Reliability  Modeling."  A.  Lemoine  and  M. 

Wenocur. 

4.  "Brownian  Motion  with  Quadratic  Killing  and  Some  Implications," 

M.  Wenocur. 

5.  "On  Sojourn  Time  in  Jackson  Networks  of  Queues,"  A.  Lemoine. 

6.  'iDiffusion  First  Passage  Times;  Approximations  and  Related 
Differential  Exjuations,"  M.  Wenocur. 

7.  'The  Loop  Elimination  Algorithm  for  Inveri/ing  Large  Sparse 
Matrices,"  R.  Barack  and  R.  Emberson. 

8.  'Approximating  Probability  Densities  on  the  Positive  Half-line,"  M. 

Wenocur. 

9.  "Waiting  Time  and  Workload  in  Queues  with  Periodic  Poisson  for 

Input,"  A.  Lemoine.  3FAJcl 

1.  .  .  . 

10.  A  Reliability  Model  Based  on  the  Gamma  Process  and  its  Analytic  ,uri«*d 
Theory,"  M.  Wenocur.  o  p  1 1  otL 


The  papers  listed  above  are  now  summarized  to  provide  background  and 
perspective  on  current  problems  of  interest  and  further  areas  for  investigation. 
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.Avail  aad/or 
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2.1  SUMMARY  OF  RESEARCH  REPORTS 

The  first  paper  presents  a  very  promising  approach  to  failure  modeling  which 
takes  account  of  the  dynamic  dependency  of  system  failure  and  decay  on  the  state  of 
the  system.  Under  this  approach,  system  state  or  wear  and  tear  is  modeled  by  an 
appropriately  chosen  random  process,  eg,  a  diffusion  process;  and  the  occurrences  of 
fatal  shocks  are  modeled  by  a  Poisson  process  whose  rate  function  is  state-dependent. 
The  system  is  said  to  fail  when  either  wear  and  tear  accumulates  beyond  an 
acceptable  or  safe  level  or  when  a  fatal  shock  occurs.  The  approach  has  significant 
merit.  First,  it  provides  new  and  revealing  insights  into  most  of  the  famous  and 
frequently  used  lifetime  distributions  in  reliability  theory,  including  the  Makeham, 
Gompertz,  Weibull,  Rayleigh  and  Gumbel  distributions.  In  fact,  these  classic  models 
are  obtained  in  a  unified  and  surprisingly  straightforward  manner.  More 
significantly,  however,  the  approach  suggests  intuitively  appealing  and 
computationally  tractable  ways  of  enhancing  these  standard  failure  models  and  for 
developing  new  ones.  [Naval  Research  Logistics  Quarterly,  Vol.  S2,  1985,  pp.  4^7- 
508.) 

Stochastic  networks,  or  vector  random  processes,  provide  an  appealing  theoretical 
framework  for  modeling  complex  sequential  processes  evolving  under  uncertainty. 
The  second  paper  describes  a  network  formulation  for  a  process  wherein  an  object  or 
"system"  moves  through  a  succession  of  states  (nodes)  and  operating  modes  (classes) 
in  the  course  of  carrying  out  its  function  (fulfilling  its  purpose).  Transitions  between 
states  and  operating  modes  occur  in  a  possibly  random  manner,  and  require 
(consume)  some  resource  in  randomly  varying  amounts.  The  paper  discusses  the 
routing  behavior  and  resource  requirements  of  a  typical  object  as  it  moves  through 
(and  eventually  out  of)  the  network.  We  then  shift  our  focus,  from  a  single  object 
and  its  odyssey,  to  the  network  as  whole,  where  time  is  the  resource  and  many 
objects  are  entering  the  network  according  to  a  possibly  nonhomogeneous  Poisson 
pattern;  in  this  vein  we  discuss  the  evolution  of  the  network  over  time.  Finally,  we 
consider  the  richness  of  applications  for  the  formulation  and  results.  [Naval  Research 
Logistics  Quarterly,  Vol.  SS,  1986,  pp.  4S1-44S.) 

The  third  paper  is  motivated  by  the  general  approach  to  .  lity  modeling  of 
the  first  paper,  and  explores  the  implications  of  shot-noise  foru.  .tions  for  s>'stem 
stress.  Suppose  the  system  is  subjected  to  shocks  or  jolts  according  to  a  Poisson 
process  with  rate  X.  Suppose  that  if  a  jolt  of  magnitude  D  occurs  at  epoch  5  then 
at  time  S+t  the  contribution  of  the  jolt  to  the  system  stress  is  Dk[t),  where  h[t)  is  a 
non-negative  function,  vanishing  on  the  negative  half-line.  If  {  5„,  n>l  }  are  the 
epochs  of  shot  occurrences  and  {  n>l  }  are  the  magnitudes  of  the  successive 
jolts,  the  "residual  system  stress"  at  time  t,  say  X[t),  is  given  by 

m  =  LD,h[t-s,). 

n—l 

Exponential  decay  of  individual  shocks  corresponds  to  h[t)  =  exp[—at).  The  time  to 
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system  failure  is  taken  as  the  epoch  of  the  first  count  in  a  doubly  stochastic  Poisson 
process  with  rate  function  X(t)-  This  formulation  is  intuitively  appealing  and  gives 
tractable  results.  Indeed,  a  class  of  failure  distributions  is  derived  using  the  shot  noise 
model  for  system  stress.  The  properties  of  these  distributions  and  their  structure  are 
explored.  [Operations  Research,  Vol.  $4,  1986,  pp.  S20-S2S.) 


The  fourth  paper  carries  forward  the  general  approach  to  reliability  modeling  in 
the  first  paper.  System  wear  and  tear  is  modeled  as  a  Brownian  motion  with  drift. 
Failure  is  due  to  shock  only  (the  system  dies  in  the  line  of  duty).  In  state  x,  shocks 
occur  at  rate  Xx^.  Thus,  the  time  to  system  failure  is  the  "killing  time,"  say  r,  of 
Brownian  motion  with  drift  and  quadratic  killing  rate.  The  paper  derives  an  explicit 
formula  for  the  state-dependent  survival  function  of  this  process,  that 
is,  P^{  T  >  t  }.  The  derivation  of  this  formula  is  truly  stunning,  involving  the 
Karhunen-Loeve  expansion  for  Brownian  motion,  special  function  theory,  and  the 
calculus  of  residues.  Some  properties  of  this  killing  time  distribution  are  also 
analyzed.  [Journal  of  Applied  Probability,  Vol.  23,  1986,  pp.  893-903.) 


The  fifth  paper  is  concerned  with  representations  for  equilibrium  sojourn  time 
distributions  in  Jackson  networks  of  queues.  For  a  network  with  N  single-server 
nodes  let  A,-  be  the  Laplace  transform  of  the  residual  system  sojourn  time  for  a 
customer  "arriving"  to  node  i  ,  "arrival"  meaning  external  input  or  internal  transfer. 
The  transforms  {  A,  :  i=l,...,N  }  are  shown  to  satisfy  a  system  of  equations  we  call 
the  network  flow  equations.  These  equations  also  lead  to  a  general  recursive 
representation  for  the  higher  moments  of  the  sojourn  time  variables  {7,  :  j=l,..A^}. 
Combining  these  formulas  with  an  appropriate  martingale,  explicit  expressions  are 
obtained  for  the  second  moment  of  sojourn  time,  and  the  covariance  of  sojourn  time 
and  queue  size  found  upon  "arrival,"  in  single-server  Markovian  queues  with 
feedback;  the  expression  for  second  moment  of  sojourn  time  derived  here  coincides 
with  the  result  reported  by  Takacs  in  1963.  The  paper  employs  probabilistic 
methods  to  obtain  probabilistic  results.  These  results  and  the  methods  employed 
hold  promise  for  making  further  b'-eakthroughs  on  an  Important  problem  which  has 
remained  remarkably  resistant  to  analysis.  [Journal  of  Applied  Probability,  Vol.  24, 
1987,  pp.  495-510.) 

The  sixth  paper  is  concerned  with  practical  and  theoretical  considerations  in 
computing  first  passage  time  statistics.  A  finite  spectral  expansion  employing  the 
method  of  finite  moments  is  derived  for  the  first  passage  distribution  and  is  shown  to 
converge  at  nearly  factorial  rates.  The  paper  is  motivated  by  first  passage  times  as 
models  of  failure  times.  Indeed,  this  paper  continues  the  line  of  development 
initiated  in  the  first  paper,  where  a  stochastic  pro^'ess  is  used  to  model  system  state, 
ie,  wear-and-tear,  and  the  system  fails  when  either  a  traumatic  killing  event  occurs 
(killing  events  happen  with  rate  k[x)  in  state  x),  or  the  system  is  retired  when  wear- 
and-tear  reaches  some  predefined  threshold  (ic,  a  first  passage  occurs).  For  example, 
if  syc’""-!  state  is  modeled  as  Brownian  motion  with  positive  drift,  then  first  passage 
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to  a  specified  threshold  has  an  inverse  Gaussian  distribution.  This  first  passage 
distribution  has  been  successfully  applied  to  numerous  problems  to  obtain  good  fits. 
A  related  but  parallel  line  of  development  is  explored  in  the  fourth  paper,  where  the 
killing  time  distribution  of  Brownian  motion  with  quadratic  killing  rate  is  calculated. 
The  aim  of  this  sixth  paper  is  to  study  first  passage  time  distributions,  where  the 
system  state  process  is  a  general  diffusion  with  reflection  at  the  origin  and  absorption 
at  r  <■  oa  That  is,  the  system  state  evolves  as  a  diffusion,  and  failure  occurs  at  the 
epoch  of  first  passage  (or  absorption)  to  level  r  .  {.Stochastic  Processes  and  Their 
Applications,  Vol.  27,  1988,  pp.  159-177.) 

The  seventh  paper  presents  a  promising  method  for  inverting  large  sparse 
matrices  of  form  I—P,  where  I  is  the  identity  matrix  and  P  is  nonnegative.  For  P 
sub-stochastic,  matrices  of  this  form  are  basic  to  Markovian  queueing  networks, 
which  represent  a  fundamental  approach  to  modeling  and  performance  evaluation  of 
large,  distributed  computer/communications  systems.  The  inversion  technique  of  this 
paper,  called  the  Loop  Elimination  Algorithm,  has  been  developed  and  tested 
extensively.  The  Loop  Elimination  Algorithm  is  quite  general,  however,  and  should 
find  application  to  flow  problems  in  large  networks.  (To  be  submitted  for 
publication.) 

The  eighth  paper  presents  and  discusses  a  variety  of  methods  for  approximating 
probability  density  functions  on  the  positive  half-line.  The  paper  focuses  on  the  class 
of  approximants  having  rational  Laplace  transforms.  It  is  well  known  that  such 
approximation  is  possible,  but  explicit  methods  for  obtaining  approximations  have  not 
been  extensively  studied.  A  survey  of  appropriate  numerical  analysis  techniques  is 
provided,  with  some  new  twists  relevant  to  the  goal  of  approximating  probability 
densities.  In  particular,  the  method  of  moments  and  orthogonal  expansion  methods 
are  studied.  A  new  proof  is  given  that  continuous  probability  densities  vanishing  at 
oo  can  be  uniformly  approximated  by  generalized  hyper-exponential  densities.  The 
same  denseness  property  is  also  shown  to  hold  for  families  of  densities  expressible  as 
sums  of  Erlang  densities  with  common  fixed  rate  parameter.  (Submitted  for 
publication.) 

The  ninth  paper  carries  forward  the  line  of  development  for  queues  with  periodic 
Poisson  input  initiated  in  Harrison  and  Lemoine  [d]  and  continued  in  Lemoine  [10|. 
This  paper  signals  a  breakthrough  on  a  difficult  problem  in  the  field,  which  has 
attracted  considerable  interest.  Formulas  are  obtained  for  asymptotic  workload  and 
waiting  time  distributions.  These  formulas  Involve  the  corresponding  moments  of 
waiting-time  (workload)  for  the  M/G/l  system  with  the  same  average  arrival  rate 
and  service  distribution.  In  certain  cases,  all  the  terms  in  the  formulas  can  be 
computed  exactly,  including  moments  of  workload  at  each  'time  of  day.’  The 
approach  makes  use  of  an  asymptotic  version  of  the  Takacs  integro-differential 
equation  for  the  time-dependent  distribution  of  workload  state,  together  with  results 
from  [4]  and  [lO].  (To  appear  in  Journal  of  Applied  Probability.) 
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The  tenth  paper  studies  reliability  models  using  the  gamma  process  as  a  stress 
driver,  continuing  a  line  of  research  first  suggested  in  Lernoine  and  VV^enocur  [ll], 
[12].  A  reliability  model  is  developed  where  system  state  is  a  random  process 
satisfying  a  stochastic  differential  equation  in  which  the  driving  process  is  gamma 
distributed.  The  necessary  Markov  process  analytic  theory  is  developed.  In  particular, 
the  infinitesimal  generator  for  the  state  ^^rocess  is  obtained,  and  an  integiu-diiferential 
equation  for  Feynman-Kac  functionals  Is  derived  and  shown  to  uniquely  determine 
the  functional.  (Submitted  for  publication.) 
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3.  DIRECTIONS  FOR  FURTHER  RESEARCH  ON 
PROBABILISTIC  METHODS 

Directions  for  further  research  on  probabilistic  methods  prompted  by  our  results 
and  suggested  by  preliminary  exploration,  include  the  following: 

3.1  PARAMETER  ESTIMATION  AND  DATA  ANALYSIS 

The  state-dependent  approach  of  Lemoine  and  Wenocur  [ll]  to  modeling  system 
failures  presents  both  significant  challenges  and  opportunities  for  statistical  analysis. 
In  order  to  apply  this  approach,  one  must  estimate  both  state-process  parameters  and 
the  killing  rate  function.  In  most  reliability  modeling  efforts  there  will  be  a  natural 
candidate  for  system  state,  and  the  corresponding  state  process  will  be  observable. 
The  observability  permits  the  study  of  conditional  probability  of  failure  given  the 
state  process.  Indeed,  the  conditional  distribution  of  failure  due  to  trauma  is  simply 
the  distribution  of  the  first  event  in  a  nonhomogeneous  Poisson  process  whose  rate  at 
time  t  is  equal  to  k(X(t)),  where  X(  )  is  the  system  state  process  and  k(  )  is  the 
killing  rate  function.  It  is  this  fact  that  allows  the  analyst  to  decompose  the 
parameter  estimation  problem  into  two  distinct  and  loosely  coupled  statistical 
estimation  problems;  namely,  the  estimation  of  the  killing  rate  function  and  the 
estimation  of  the  system  state  process  parameters. 

Furthermore,  with  the  evolution  of  micro-electronics  hardware,  there  is  an  ever 
increasing  system  self-monitoring  capability.  The  benefits  of  self-monitoring  can  only 
be  properly  exploited  by  using  state  dependent  failure  modeling.  A  great  deal  can  be 
accomplished  using  data  collected  by  self-monitoring  systems.  It  can  be  used  to 
improve  understanding  of  failure  as  a  function  of  system  state.  Such  understanding 
can  result  in  better  methods  for  early  failure  prediction.  An  important  first  step  in 
developing  the  mathematical  techniques  needed  to  exploit  monitoring  information  is 
embodied  by  the  dynamic  failure  modeling  approach  suggested  in  [ll]. 

We  plan  to  explore  estimation  problems  suggested  by  the  state-dependent 
approach  to  mo  (eling  system  failures,  that  is,  estimation  of  state-process  parameters 
and  killing-rate  function.  A  starting  point  is  the  literature  on  statistical  inference  for 
stochastic  processes,  cf.  Basawa  and  Prakasa  Rao  [2],  Kutoyants  [9],  and  Nanthi  [18]. 
We  anticipate  that  methods  of  conditional  maximum  likelihood,  cf.  Cox  and  Hinkley 
[3];  survival  analysis,  cf.  .Miller  [17];  multiplicative  intensity  models  for  multivariate 
counting  processes,  cf.  Jacobsen  [5\;  point  processes,  cf.  Karr  [7j;  and  martingales 
and  stochastic  integrals,  -'f.  Metivier  and  Pellaumail  [15],  Metivier  [16|,  and  Karatzas 
and  Shreve  [6],  will  be  useful. 


3.2  SOJOURN  TIME  IN  NETWORKS  OF  QUEUES 

Based  on  the  results  of  Lemoine  [13],  we  envision  several  promising  approaches  to 
approximating  the  higher  moments  of  sojourn  times  in  Jackson  networks  of  queues. 
In  particular  we  have  derived  a  functional  equation  which  permits  us  to  approximate 
the  second  moments.  This  functional  equation  has  some  intriguing  properties,  eg,  it  is 
a  positive  linear  map  from  /[(tr)  sequences  to  /j(7r)  ,  indeed  a  contraction  mapping 
possessing  a  unique  solution.  We  are  cautiously  optimistic  that  this  approach  could 
lead  to  a  closed  form  solution  for  the  second  moments  of  sojourn  times. 

We  believe  also  that  martingale  methods  and  time  reversal  combined  with  our 
representation  results  for  stochastic  networks  could  provide  the  means  to  obtain 
approximations  (or  exact  solutions)  to  higher  moment  problems  and  Laplace 
transforms. 

3.3  MOMENT  CONSTRAINED  ADAPTIVE  METHODS 
FOR  ORTHOGONAL  EXPANSIONS 

This  method,  suggested  in  Section  5  of  Wenocur  [21],  Holds  great  promise  for 
improving  the  accuracy  of  probability  density  approximation  on  the  positive  half-line. 
Before  this  method  can  produce  useful  results,  a  number  of  important  issues  must  be 
investigated.  Namely,  find  conditions  under  which  the  supremum  is  obtainable  and 
unique;  identify  what  special  structure  may  be  exploited  to  obtain  efficient 
algorithms;  and  how  stable  and  computationally  expensive  this  approach  may  be. 

We  envision  this  task  will  require  both  theoretical  and  computational  efforts. 


3.4  IMPROVED  BERNSTEIN  TYPE  APPROXIMATIONS 

Bernstein  polynomial  approximation  is  an  important  method  for  both  theoretical 
and  numerical  analysis.  Bernstein  polynomial  approximation  being  a  positive 
operator  makes  it  very  attractive  as  a  technique  for  probability  density 
approximation  (cf.  Section  2  of  VVTnocur  [21]).  However,  its  slow  rate  of  convergence 
is  a  significant  drawback  in  this  context.  A  theorem  of  Korovkin  (cf.  Chapter  -1  of 
[8])  states  that  positive  linear  operators  can  achieve  n~^  rate  of  convergence  but  no 
more.  In  contrast,  standard  Bernstein  polynomials  achieve  only  rate  of 

convergence  . 

Our  goal  is  to  produce  an  approximation  method  similar  to  Bernstein 
approximation,  but  obtaining  the  optimal  n'"  optimal  rate.  Our  preliminary  study 
has  produced  a  method  having  n~'^^'Uog{n)  rate  of  convergence. 
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3.5  IMPROV]  DIFFUSION  APPROXIMATIONS  FOR 
BIRTH-DUaTH  MARKOV  CHAINS 

Preliminary  results  strongly  suggest  that  alternate  parametrization  may 
significantly  improve  accuracy  of  diffusion  approximations  for  M/M/s  under  light  f)r 
heavy  traffic  conditions  when  compared  to  the  standard  approach. 


3.6  RELIABILITY  MODELS  BASED  ON  DIFFUSION  AND 
OTHER  MARKOV  PROCESSES 

We  intend  to  carry  forw'ard  the  line  of  development  for  dynamic  reliability 
modeling  initiated  in  Lemoine  and  Wenocur  [ll]  and  continued  in  Wenocur  [19],  ['i'ij. 
In  particular,  we  will  extend  the  results  of  Wenoctir  [19]  to  the  setting  where  the 
Ornstein-Uhlenbeck  process  is  the  model  for  system  stress  and  traumatic  killing 
events  are  a  quadratic  function  of  system  state.  We  will  also  investigate  theoretical 
and  computational  aspects  of  the  Markov  analytic  theory  developed  in  [22|. 


3.7  QUEUES  WITH  PERIODIC  POISSON  INPUT 

As  shown  in  [14],  exact  computation  of  moments  for  asymptotic  workload  and 
waiting  time  distributions  involves  the  values  {/4(0)}  where  //,(0)  is  the  long-run 
probability  that  server  workload  is  zero  at  ’time  of  day  s.'  When  these  probabilities 
are  known  the  moments  can  be  calculated  explicitly.  If  the  service  distribution  is 
concentrated  on  integer  multiples  of  the  period  length,  these  probabilities  are  equal 
(over  the  period)  and  known.  Our  objective  is  to  compute  (or  ap[)roximate)  the 
values  {Hg{0)}  for  arbitrary  distributions  of  service  time,  and  to  explore  implications 
for  calculating  (or  estimating)  asymptotic  distributions  for  workload  and  waiting 
time. 


3.8  INVERTING  LARGE  SPARSE  MATRICES 

Our  goal  is  to  characterize  the  class  of  matrices  for  which  the  L(X)p  I'dimination 
.‘Mgorithm  [l]  converges,  and/or  find  an  efficient  method  for  checking  same.  (.Vn 
appropriate  analogy  is  cycling  in  the  simplex  method.) 
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4.  RELATED  CONTRACT  ACTIVITY 

Supported  in  part  by  the  subject  AFSC  contract,  A.  Lemoine  and  M.  Wenocur 
attended  the  following  professional  meeting; 

•  Special  Interest  Conference,  "Queueing  Networks  and  Their  Applications."  This 
conference  was  held  at  the  Hyatt  Regency  Hotel  in  New  Brunswick,  New  Jersey, 
January  7-9,  1987. 

•  Sixteenth  (16th)  Conference  on  Stochastic  Processes  and  their  Applications.  This 
conference  was  held  at  Stanford  University.  .-Vugust  16-21,  1987. 
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